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ABSTRACT 

Solar flares involve complex processes that are coupled and span a wide range of temporal, spatial, and en¬ 
ergy scales. Modeling such processes self-consistently has been a challenge in the past. Here we present results 
from simulations that couple particle kinetics with hydrodynamics of the atmospheric plasma. We combine the 
Stanford unifled Eokker-Planck code that models particle acceleration and transport with the RADYN hydro- 
dynamic code that models the atmospheric response to collisional heating by accelerated electrons through 
detailed radiative transfer calculations. We perform simulations using two different electron spectra, one an ad 
hoc power law and the other predicted by the model of stochastic acceleration by turbulence or plasma waves. 
Surprisingly, the later model, even with energy flux 10^^ erg s“^ cm“^, can cause “explosive” chromospheric 
evaporation and drive stronger up- and downflows (and hydrodynamic shocks). This is partly because our 
acceleration model, like many others, produces a spectrum consisting of a quasi-thermal component plus a 
power-law tail. We synthesize emission line profiles covering different heights in the lower atmosphere, in¬ 
cluding Hof 6563 A, He II 304 A, Ca II K 3934 A and Si IV 1393 A. One interesting result is the unusual high 
temperature (up to a few 10^ K) of the formation site of He II 304 A, which is expected due to photoionization- 
recombination under flare conditions, compared to those in the quiet Sun dominated by collisional excitation. 
When compared with observations, our results can constrain the properties of non-thermal electrons and thus 
the poorly understood particle acceleration mechanism. 

Subject headings: acceleration of particles — hydrodynamics (HD) — line: profiles — radiative transfer — 

Sun: flares; chromosphere 


1. INTRODUCTION 

One of the outstanding problems in solar physics is how 
magnetic energy is transformed into the observed signatures 
of solar flares. It has been recently recognized that model¬ 
ing the coupling between the particle acceleration and trans¬ 
port processes in solar flares and the dynamical response of 
the atmosphere to particle collisional heating is critical to 
our understanding of flare dynamics. It is clear that non- 
thermal electrons and ions play an important role. How¬ 
ever, the exact mechanism of acceleration of these particles 
is still a matter of considerable debate. Several scenarios 
have been proposed and different models have been devel¬ 
oped with different degrees of details . Among thes e are ac- 
celer ation by DC electric fie lds (e.g.. lHolmanllI985l) . shocks 
(e.g.. iTsuneta & NaitalI998h . and turbulence. Stochastic ac¬ 
celeration by turbulence or plasma wave s (iRamatvl II979I: 
iHamilton & PetrosianI r 9921 iMilleil Il997h ha s been devel¬ 
oped in greater detail (iPetrosian & Lid I2Q04I) and has been 
tested by observation s more rigorously (IPetrosian & Donaghv 
119991: [Liu et al .11200^ l2006b : for a recent review see IPetrosian 

Observations most intimately connected to the acceleration 
process are the microwave. Hard X-ray (HXR) and y-ray radi¬ 
ations. These emissions are produced by particles, which are 
(most likely) accelerated in the corona in a flaring loop. The 
accelerated electrons produce microwave and HXR emissions 
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via synchrotmn and bremsstrahlung p rocesses (iHovng et al.l 
II98u IS^^II994l: iKundu et al.llI994)) . while the interaction 
of accelerate d protons a nd ions with the background ions pro¬ 
duce y-rays (lLinllI985l) . However, most of the energy of the 
particles goes into heating of the plasma via Coulomb col¬ 
lisions as they travel down to the footpoints. The heating 
and evaporation of the plasma produces secondary continuum 
and line emissions, from infrared to soft X-ray energy bands, 
which also carry information about the acceleration process. 
The purpose of the work presented here is to explore this sec¬ 
ond avenue of testing the acceleration models. 

Eor the long run, our aim is to use the combined non- 
thermal and thermal signatures to distinguish between the 
different acceleration scenarios and constrain the character¬ 
istics of the specific acceleration models. This requires a 
combined treatment of the acceleration, transport and radi¬ 
ation of particles and hydrodynamic (HD) response of the 
atmosphere to the energy input by particles. The first nu¬ 
merical study_oUparticle transport in solar flares was car¬ 
ried by iLeach & PetrosianI (II98II) who treated the transport 
and radiation by electrons using the Eokker-Planck transport 
equation taking into account the pitch angle changes due to 
Coulomb collisi ons and magnetic held variation s. This study 
was extended bv lMcTiernan & PetrosianI (II99QI) . who consid¬ 
ered energy loss and pitch angle changes due to synchrotron 
emission. This transport-radiation c ode was later combined 
with an stochastic acceleration code (jPark & PetrosianllI995l: 
IPetrosian & Donagh^lI999l: IPetrosian & Liul I2004I) into the 
unifled Stanford code (IPetrosian et al.ll2002lk 

The HD response has been investigated in several works 
(iMariska et al.lll989l: iKasparova et aklUOP^ by the means of 
ID numerical hydrodynamics along a coronal loop. These 
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works assumed that non-thermal electrons with a power-law 
spectrum were injected at the apex of the loop and used to 
approximate analytic expressi ons for the transport and en- 
ergy deposit i on alo n^ the loop. lAbbett & HawlevI (119991) and 
I Allred et ^ (l2005h improved the results by including a de¬ 
tailed calculation of radiative transfer in the atmosphere. 

We have embarked on development of more complete and 
self-consistent treatment of this important problem. Instead 
of using an ad-hoc power-law injected electrons and approx¬ 
imate treatment of transport, we have combined the Stanford 
Fokker-Planck acceleration-transport code with HD codes, 
achieving a more accurate determ ination of the rad iative sig¬ 
natures of flares. In our first paper. lOu et al.l (120091) (hereafter 
paper I), we combined the Stanford co de with the Naval Re - 
search Laboratory (NRL) hydro-code (iMariska et al.lll98^ . 
where we addressed some of the non-thermal aspects of the 
problem. Here we extend this study by taking into account 
the detailed calculation of radiative transfer in the atmo¬ 
sphere. Instead of the NRL code we use the radiative HD 
(RADY N') code dCarlsson & SteinlfT^fTW^ . in a modified 
version (lAbbett & HawlevI 1 19991: lAllred et al.ll2005h . and fo¬ 
cus on intensities and shapes of several lines emitted at dif¬ 
ferent heights. The main goal here is to show the effects of 
more realistic accelerated spectrum of electrons as compared 
to those of power-law injected electrons. In a recent paper 
we have used the RADYN code to interpret spectroscopic 
observations of IR-optical-U V continuum and line emissions 
(iRubio da Costa et 1201 51) . More such comparison with 
specific observations, both line and (thermal and non-thermal) 
continuum emissions will be presented in future publications 
using this combined code. 

This paper is organized as follows. In Section [2] we de¬ 
scribe our model and the method used to solve the combined 
radiative HD and Fokker-Planck equations. In Section [S] we 
present results on plasma characteristics from the simulations 
and compare them with those from the power-law injection 
case. In Section |4] we focus on the effects of electron heat¬ 
ing on the emission in the Ha, Ca II K, He II 304 A and 
Si IV 1393 A lines. A brief summary and conclusions are 
given in Section[5l In Appendix^we present details of differ¬ 
ent energy terms and their influence on the atmospheric evolu¬ 
tion and in Appendix iBl we investigate the effects of different 
electron fluxes. 

2. MODEL AND METHOD 
2.1. Assumptions and Geometry 

Magnetized plasma processes, in general, are com¬ 
monly modeled with multidimensional magnetohydrodynam¬ 
ics (MHD) codes. However, in the solar corona the magnetic 
pressure dominates over the gas pressure (i.e. the plasma 
P 1), so plasma flows mainly along magnetic field lines. 
This allows us to use the HD instead of MHD treatment of the 
problem. 

We assume a one-dimensional, semi-circular loop perpen¬ 
dicular to the solar surface with a constant diameter of 3 Mm. 
Adopting a symmetric boundary condition at the loop apex, 
we construct a computational domain containing a 10.4 Mm 
long quarter circle, discretized in 191 grid points (a denser 
grid would be computationally too expensive). The loop ex¬ 
tends in a plane-parallel model atmosphere from the corona to 
the bottom of the chromosphere at z = 0 Mm, where the op¬ 
tical depth T 5000 at a wavelength of T = 5000 A is unity. The 
circular geometry is taken into account when calculating the 


X-ray photoionization rates and the gravitational acceleration; 
otherwise the loop is treated as a vertical cylinder. 

For the pre-fla r e con ditions we use the FP2 model of 
lAbbett & HawlevI (119991) . which was generated by adding 
a transition region and corona to the model atmosphere of 
ICarlsson & SteinI (Il997h . The temperature is fixed to 10^ K 
at the loop apex. By running the code without external heat¬ 
ing, the atmosphere relaxes to a state of HD equilibrium. 

2.2. Radiative Transfer and Hydrodynamics 

The RADYN code simultaneously solves the equations of 
hydrodynamics, population conservation, and mdiative trans¬ 
fer u sing a one-dimensional adaptive grid (iDorfl & DrurvI 

Atoms important to the chromospheric energy balance are 
treated in non-local thermodynamical equilibrium (non-LTE). 
These include a six-level plus continuum hydrogen atom; 
a six-level plus continuum, singly ionized calcium atom; 
a nine-level plus continuum helium atom; and a four-level 
plus continuum, singly ionized magnesium atom. The tran- 
sitions that are treated in detail are given in Table 1 of 
lAbbett & HawlevI (119991) . Complete redistribution is assumed 
for all lines, except for the Lyman series in which partial fre¬ 
quency re-distributio n is mimicked by truncati ng the profiles 
at 10 Doppler widths (iMilkev & Miha^lasll 19731) . Other atomic 
species are included in the calculation as backgrou nd continua 
in LTE , using the Uppsala opacity package of iGustafssonI 

(CItI. 
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Figure 1. Spectra of angle-integrated electron fluxes, F{E), times at the 
loop-top and tmax-^^ s for different simulation runs labeled with the corre¬ 
sponding energy fluxes The black line represents a power law with 

a spectral index d = 5 and low energy cutoff Ec = 15 keV (Run PL), 
while the colored lines represent stochastically accelerated electron spectra 
(Runs SA1-SA3). [A color version of this figure is available in the o nlin e 
journal.]. 

As described in lAllred et al.l (l2005h . the radiative HD code 
includes the calculation of photoionization heating resulting 
from high temperature, soft X-ray emitting regions, as well 
as the calculation of optically-thin cooling due to thermal 
bremsstrahlung and collisionally excited metal transitions. 
An adjustment in the calculation of the conductive flux has 
also been taken into account in order to avoid unphysical large 
values in the transition region - where temperature gradients 
are large. Hydrodynamic effects due to gravity, thermal con- 
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Table 1 

Summary of the different simulation runs and the atmospheric parameters for different flaring conditions. 


Run 

Injected e distribution 

Smax 

(erg s“^ cm“^) 

^max 

(km s“^) 

^(^max) 

(s) 

^min 

(km s“^) 

K^min) 

(s) 

Tmax 

(K) 

K^max) 

(s) 

PL 

d=5;£c=15keV 

1.2x10“' 

483 

46 

-44 

30 

9.5 X 10® 

75 

SAl 

Stochastic acceleration 

1.2x10'® 

750 

32 

-41 

10 

2.3 X lO'' 

60 

SA2 

Stochastic acceleration 

6.2x10® 

641 

37 

-39 

13 

1.8 X 10’ 

60 

SA3 

Stochastic acceleration 

5.7x10* 

377 

77 

-28 

31 

8.0 X 10® 

61 


Notes: i^max - maximum electron energy flux. Vmax & K^max)^ ^min & and T^^x & KTmax) are the maximum velocity (upflow, v > 0), minimum velocity 

(downflow, V < 0), and maximum temperature and their corresponding times, respectively. 


t = 5 secs t = 20 secs t = 60 secs t=120 secs 



Figure 2. Variation with distance along the loop of the atmospheric parameters (temperature, electron density, velocity of the plasma and electron heating rate 
per unit volume (Qe) where positive velocity refer to upflows and negative velocities to downflows) in the PL (red line) and SAl (blue line thick) models. The 
black dashed line show the initial values and the dotted lines in the bottom panels represent the sound speed. The columns are for r = 5, 20, 60 and 120 s. (see 
the o nlin e movie for more details.). [A color version of this figure is available in the online journal.]. 


duction, and compressional viscos ity are considered as de¬ 
scribed in lAbbett & Hawle^^ (119991) . 

Inclusion of the radiative transfer calculation gives us the 
advantage, with respect our Paper I, to investigate how the 
electron deposition in the chromosphere affects to the emis¬ 
sion lines originating from several heights (see Section lU). 

2.3. Electron Acceleration, Transport and Energy Deposition 

FLARE - the Stanford unifi e d acceleration-transpor t code 
(iMcTiernan & PetrosianI 1 19901: iPetrosian & Liul 120041) con¬ 
sists of two modules. The acceleration module calculates the 
spectrum of the electrons accelerated stochastically by tur¬ 
bulence in the acceleration region (assumed to be located at 


the lo op-top; for observational evidence see, e.g., iLiu et al.l 
12013h and the spectrum of electrons escaping down to the 
footpoints. The transport module then calculates the evolution 
of the spectrum and pitch angle distribution along the loop of 
the escaping electrons from the acceleration site. The code in¬ 
cludes energy and pitch-angle diffusion due to Coulomb col¬ 
lisions and synchrotron radiation. 

One of the main inputs to the RADYN code is the heat¬ 
ing rate Qe{s) (in units of erg s“^ cm“^) due to acceler¬ 
ated electrons as a function of position 5* along the loop, 
which is included as a source of extern al he ating in the en- 
ergy conservation equatio n (see Equation lATl in AppendixiAb. 
lAbbett & HawlevI (119991) assumed a beam of electrons with 
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a power-law spectrum with an index S and low energ y cutoff 
Er and used the approximate analytic expression of lEmslid 
(1198 Ih to calculate the heating rate. Here, we use the electron 
spectrum from the acceleration-transport code, which in gen¬ 
eral includes a quasi-thermal plus a non-thermal component. 
We then calculate the heating rate Qe{s) using the spatial vari¬ 
ation of the spectrum and pitch angle distribution f(E,ju, s) 
of the electrons, following the procedure described in Sec¬ 
tion 3.1 of Paper I. 

2.4. Combining ELARE and RADYN 

We followed the approach detailed in Paper I to combine 
the FLARE and RADYN codes. In brief, we note that the 
current particle transport module of FLARE can provide only 
a steady-state solution. To perform time-dependent simula¬ 
tion, RADYN calls the transport module every At = 1 s. This 
interval of communication between the two codes allows for a 
steady-state solution of particle kinetics, whose timescales are 
orders of magnitude shorter than HD timescales. This interval 
also provides a trade-olf for computational efficiency. 

In each iteration, the RADYN code passes the electron den¬ 
sity of the atmosphere as a function of distance from the top of 
the loop (ne(s)) to the FLARE code. The FLARE code calcu¬ 
lates the electron heating rate Qe(N) as a function of column 
density Ne = ne(s) ds and return it to the RADYN code, 
which then simulates the evolution of the atmosphere and ad¬ 
vance in time for Ar = 1 s. Within this interval, we assume 
that the accelerated electron flux remains constant and so does 
QeiNe). Therefore, the temporal evolution of the heating rate 
Qe{s, t) as a function of distance solely depends on the spa¬ 
tial redistribution of the electron density ne(s, t) with time, as 
described in Section 3.2 of Paper I. 

By including the radiative transfer calculation, one of the 
advantages with respect our previous Paper I is the study in 
detail how the electron deposition in the chromosphere affects 
to the emission of the lines at several heights (see Section lU). 

2.5. Simulation Runs 

We have performed four simulation runs with different in¬ 
jected electron spectra, as summarized in Table [TJ to inves¬ 
tigate their different atmospheric response. These runs were 
carefully designed and allowed us to compare our results us¬ 
ing the acceleration model based spectra with those carried 
out in the past using a simple power-law. In the Run PL, we 
use a single power-law electron spectrum of index 6=5 and 
lo w energy cutoff of £^ ^=15 keV and the analytic expression 
of lEmslid (ll978L[T9Mh to calculate the electron heating rate 
Qe. Runs SA1-SA3 use the result of the stochastic acceler¬ 
ation (SA) model for the injected electron spectra, for which 
we prescribed a common characterist ic acceleration timescale 
of Tp = 1/70 s (iPetrosian & Liull2004l) and thus the same spec¬ 
tral shape. For the SA runs we c alculate the Qe using the 
approach described in Section [231 Note that Run SAl has 
identical electron spectra with Run N in Paper I and Run PL 
has a lower spectral index than Run O in Paper I. Each simu¬ 
lation lasts 120 s with the electron energy flux g(0 increasing 
linearly till r^^_^=60 s up to gmax and then decreasing for an¬ 
other 60 s. In all the runs we keep the electron spectral shape 
constant in time. 

Figure [T] compares the spectra of angle-integrated electron 
flux at the top of the loop for the different runs shown in 
Table [T] In the spectral range approximated for hard X-ray 
bremsstrahlung emission, the slope of the PL spectrum has 


been chosen to be similar to the one of the SA Run. The 
main difference between the two models is the absence of the 
quasi-thermal electrons in the PL model and the secondary 
difference is the roll over in the SA model at highest energies. 

We should, however, note that these differences are not 
unique to the SA model. Most models, when accelerating 
particles from a thermal pool, produce a hotter quasi-thermal 
component with a power law tail because of Coulomb col- 
lisi on effects. Interested re aders are referred to, for exam¬ 
ple, Petrosian & Fasti (120081) for a generic acceleration model 
and lSironi & Spitkovskvl (120091) for shock acceleration using 
particle-in-cell (PIC) simulations. We also note that the selec¬ 
tion of the spectral parameters has not been done to study any 
particular flare. 

In the next section we investigate how the atmosphere re¬ 
sponds to the injection of single power law electrons (Run PL) 
and stochastically accelerated electrons (Run SAl) and in Ap¬ 
pendix [B] we compare results of Runs SAl, SA2 and SA3 to 
investigate how the atmosphere responds to the variation of 
the electron flux. 

3. COMPARISON OF SINGLE POWER-LAW INJECTION WITH 
STOCHASTIC ACCELERATION 

In order to investigate how the acceleration and transport 
of electrons affects the atmospheric response, we compare the 
results of injection of single power-law electrons (Run PL) 
and of stochastically accelerated electrons (Run SAl), which 
have the same total energy flux at any moment. 

3.1. Atmospheric evolution 

Figure [2| and the online movie shows the temporal evolu¬ 
tion of the atmosphere. In general SAl has higher values than 
PL for all the atmospheric variables (see Table [T]), qualita¬ 
tively consistent with those of Runs N and O in Paper I, re¬ 
spectively. For example, the corona is heated more rapidly to 
higher temperatures in SAl than that in PL. The maximum 
coronal temperature of 2.3 x 10^ K is reached at r = 60 s in 
SA, compared with a lower value of 9.5 x 10^ K at 74 s in PL 
(Figure Oc)). As detailed in Section [3^ such contrasts are 
due to the different spatial distributions of collisional heating 
by non-thermal electrons of the two models, even though they 
have the same energy input to the loop as a whole. This factor 
of 2.4 difference in the maximum coronal temperature is much 
greater than the factor of 1.2 difference between similar runs 
in Paper I (see their Table 1), which points to the importance 
of detailed radiative transfer calculation included here. 

At early times, electron heating in the transition region and 
chromosphere causes an overpressure that drives both upflows 
(i.e., chromospheric evaporation) and downflows of plasma, 
with higher speeds in the case of SAl. At r = 5 s, for ex¬ 
ample, SAl has the maximum upflow and downflow veloc¬ 
ities of 102 and -7 km s“^ respectively, compared with 13 
and -0.4 km s“^ in PL (see Figure [3] for the temporal evolu¬ 
tion). As a result of the upflow, part of the initial transition- 
region and chromospheric material is converted into coronal 
mass, and the new transition region recedes to lower heights 
between z = 0.99 and 1.10 Mm. 

At r = 5 s most of the energy is deposited in a narrow region 
located in the upper chromosphere at a height of z = 1.27 Mm 
for PL and at 1.34 Mm for SAl. The non-thermal electrons 
quickly heat this region to temperatures greater than 10^ K. 
The location of the maximum non-thermal energy deposition 
gradually moves downwards, reaching a height of 1.17 Mm 
for PL and 1.12 Mm for SAl at ^ = 20 s. A more detailed 
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Figure 3. Temporal evolution of the maximum (a) upflow velocity, (b) down¬ 
flow velocity, and (c) temperature within the loop for Runs PL (red line) and 
SAl (blue line thick). The dashed curve indicates the sound speed at the cor¬ 
responding location of the maximum upflow (a) and downflow (b) velocity 
for each run. [A color version of this flgure is available in the online journal.]. 

discussion on how the electron heating rate affects the atmo¬ 
sphere is given in Appendix El where the height distribution 
of the different energy terms and their temporal evolution are 
presented. We note that in the case of PL a significant frac¬ 
tion of hydrogen is ionized at f = 5 s and z = 1.18 Mm and 
thus the energy can no longer be effectively radiated away by 
bound-bound transitions. 

Figure [3a) shows the temporal evolution of the maximum 
velocity along the loop (solid lines) and the sound velocity 
at the same position (dashed lines). Again, SAl has higher 
plasma velocities, creating a hydrodynamic shock with super¬ 
sonic speeds between t = 22 and 49 s, while in the PL case, 
a shock is developed at later times between t = 30 and 78 s. 
For SAl, the maximum upflow velocity of Vmax = 750 km s“^ 
with a Mach number of 1.15 is attained at f = 32 s, while for 
PL the maximum velocity Vmax = 483 km s“^ with a Mach 
number of 1.41 occurs at f = 46 s. 

As shown in Figure I3b), between t = 9 and 24 s for the 
SAl Run, the downflow plasma velocity exceeds the sound 
speed and forms a downward propagating shock. In contrast, 
the sound speed (not shown) for the PL case in (b) is alway 
greater than the maximum downflow speed, except for a brief 


period of 3 s from r = 32 to r = 34 s. 

In Figure |3c) we can also see that the temporal evolution 
of the maximum temperature along the loop follows a similar 
general trend as the electron flux 5(0. but the temperature in 
the SAl run increases much faster than that in the PL run. 

Figure [4] shows the temperature distribution of the plasma 
velocity at different times. The two runs exhibit distinct be¬ 
haviors at early times. At r = 5 and 20 s (panels (a) and (b)), 
for example, in the PL case, mainly low-speed upflows in the 
10^-10^ K range are present. In the SAl case, downflows oc¬ 
cur at low temperatures and high-speed upflows at high tem¬ 
peratures, with their division temperature increasing with time 
from a few 10"^ K to nearly 10^ K. Comparably slower down¬ 
flows take place in the PL case only at later times (e.g., 60 s). 
As time progresses toward the late phase of the flare, such a 
distinction diminishes to lesser degrees, with mainly upflows 
present in both cases whose speed rapidly grow with temper¬ 
ature. 

Such distinct temperature distributions of the plasma ve¬ 
locity, especially early during a flare, are manifested in the 
Doppler shifts of emission lines shown in Figures [6l [9l [l2| 
and [141 which will be discussed later, and demonstrates their 
sensitive dependence on the spectra of accelerated electrons. 
When compared with observations, such distributions can be 
used as diagnostics to constrain particle accel eration mecha¬ 
nisms. For example. iMilligan & DennisI (120091) used the EUV 
Imaging Spectrometer (FIS) onboard Hinode to measure the 
Doppler velocities of emission lines at formation tempera¬ 
tures ranging from 0.05 to 16 MK. They found a temperature 
distribution of Doppler velocity akin to that of SAl shown in 
Figure [3b). Specifically, their -60 to -30 km s“^ downflow 
speeds within 0.6-1.5 MK are comparable to those in our SAl 
case at r = 9-12 s. They found slightly higher division tem¬ 
perature of ~2 MK between up- and downflows. This minor 
difference could be partly due to the greater electron energy 
flux of 5 X 10^^ erg cm“^ s“^ estimated in their C-class flare, 
about five times that adopted in our simulation. 

3.2. Evolution of the electron heating rate 

The height distribution of the electron heating rate Qe{s) is 
shown in the bottom row of Figure[2]at t = 5, 20, 60 and 120 s 
for PL (red line) and SAl (blue thick line). In general, both 
runs present a peak in the upper chromosphere. However, Qe 
in PL is higher in the chromosphere but lower in the corona 
than that in SAL 

In PL the electron energy is mostly deposited in the upper 
chromosphere (at z = 1.18 Mm for t = Is). The electrons 
interact with the plasma, ionizing hydrogen and helium. A 
step in the heating rate appears at the chromospheric evapo¬ 
ration front due to the density jump (e.g., at z ~ 3.7 Mm in 
Figure [3n))- During the decay phase of the flare, the heating 
rate decreases throughout the loop, but there is more fractional 
reduction in the chromosphere than in the corona. 

In the case of SAl, Qe is mostly concentrated around z 

1.45 Mm after r = 1 s, covering a broader region than in PL. 
This is because the stochastically-accelerated electrons have a 
quasi-thermal spectral component plus a non-thermal tail cov¬ 
ering a broader energy range, as shown in Figured] The quasi- 
thermal component carries considerable energy contents and 
causes significant heating in the corona that reaches temper¬ 
atures of 11 MK at r = 10 s and 23 MK at 60 s. The sec¬ 
ondary heating peak in the lower chromosphere is due to the 
non-thermal tail that exceeds the power-law spectrum from 
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(a) t = 5 sec (b) t = 20 sec (c) t = 60 sec (d) t=120 sec 



Figure 4. Temporal evolution of the velocity vs. temperature for PL (red solid line) and SAl (blue solid thick line). [A color version of this hgure is available in 
the online journal.]. 


(o) Energy Spectrum (b) Spatial Distribution 



Figure 5. Energy and spatial distributions of non-thermal electrons at tmax = 60 s in Run SAl. (a) Energy spectra of angle-integrated electron flux multiplied by 
the energy squared. The solid black line is the flux {Fesc) of electrons escaping the acceleration region at the loop apex that serves as the injection to the transport 
portion of the loop. The color-coded curves (from red to blue) of different line styles represent fluxes at growing distances of 4, 6, 7, 8, 9, and 10 Mm from the 
loop apex, (b) Spatial distribution of angle-integrated electron flux as a function of distance along the loop. The color-coded lines represent fluxes at different 
electron energies, 3, 6.1, 12.3, 24.5, 48.8, 97.4, and 294.1 keV. Note the large step at the transition region, especially at low energies. [A color version of this 
figure is available in the online journal.]. 


?^800 keV to 6 MeV. The sharp spike in the photosphere is 
due to the increase of the electron density at that height. 

The above behaviors of the electron heating rate can be 
better understood by examining the energy and spatial dis¬ 
tributions of the non-thermal electrons. Figure [3 a) shows the 
energy spectra of E^F(E, s) at selected distances 5* from the 
top of the loop for Run SAl at f = 60 s, where E(E, s) = 

V f(E, yu, s)d/j is the pitch-angle integrated electron flux. 
As expected, low-energy electrons suffer most of the losses 
at higher altitudes than the high energy electrons. This is 
mainly because of the 1 / v dependence of the Coulomb colli¬ 
sion energy-loss rate. Note that the large decrease of electron 
flux from 5' = 9 to 10 Mm results from the change of loca¬ 
tion from above to below the transition region with a sudden 


increase of the atmospheric density. 

Figure Ob) shows the spatial distribution of the electron 
flux E(E, z) along the loop at selected electron energies for 
Run SAl at f = 60 s corresponding to Figure Oa). In gen¬ 
eral, the electron flux decreases with distance from the in¬ 
jection site at z = 10.40 Mm. The slope (dE(E,z)ldz) is 
steeper at lower energies because low-energy electrons loose 
energy faster than high-energy electrons. This slope also 
depends on the ambient density Ue, because dE{E,z)ldz = 
ne[dE(E,N)/dN], where dE( E, N)ldN is generally a smoot h 
function of column density N (iMcTiernan & Petrosianll 19901) . 
As such, a sharp drop in the electron flux, more prominent at 
lower energies, occurs at the transition region (z = 0.8 Mm) 
due to the sharp ambient density rise. In addition, the slope 
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Figure 6. Ha 6563 A line profiles at r = 5, 20, 60 and 120 s for runs PL (red line) and SAl (blue thick line). The emission at the initial time has been subtracted 
and the intensity has been normalized to the continuum. The vertical dashed line indicates the line center rest wavelength. 
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Figure 7. Intensity contribution function, as dehned in equation {TJ, for Her as a function of frequency (bottom axis) and height (vertical axes), shown as reversed 
color maps (darker for higher values). The frequencies are in velocity units such that positive and negative velocities represent upfiows and downflows of plasma, 
respectively. The dashed blue line represents the atmospheric velocity stratiheation; the solid green line the height at which Ty = 1 and the black thick line shows 
the line profile (including the quiet Sun emission). Top panel: SAl Run; bottom panel: PL Run. All panels have the same color scale as indicated on the colorbar. 
[A color version of this figure is available in the online journal.]. 


is proportional to the electron heating rate (see equations (10) 
and (11) in Paper I) and thus accounts for its shape, espe¬ 
cially the narrow peak at the transition region, as shown in 
Figure [3m-p). 

4. LINE EMISSION 


One o f the advantage s of this work over previous stud¬ 
ies te.g.. iLiu et al.ll200^ is the inclusion of the detailed ra¬ 
diative transfer calculation of emission lines treated in non- 
LTE and thus the capability of synthesizing line emission 
from simulation results at different height formations with 
the aim of studying how the energy deposition affects the 
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lower atmosphere. In this section we will examine four 
emission optically thick lines, Ha, Ca II K, He II 304 A, 
Si IV 1393 A, which are common observables and cover for¬ 
mation heights from the upper photosphere to the transition 
region. For example, the Ha and Ca II K lines are routinely 
observed by ground-base d facilities, e.g., at the Big Bea r So¬ 
lar Observatory (BBSO, iJohannesson et al.iri995Ul998h and 
the Kanzelhohe Observatory (KSQ. iPotzi et a]Tl2013h : the 
Extreme-ultraviolet Variability Experime nt (EVE) on board 
the Solar Dynamics Observatory {SDO) (I Woods et al.ll2012h 
covers the He II 304 A and Si IV 1393 A lines; the latter 
is also covered by the Interfac e Region Imaging Spectrograph 
{IRIS) (iDe Pontieu et al.ll2Qa . We intend to perform detailed 
comparison with observations in future investigations. 

In order to better understand the behavior of these lines, we 
write the formal solution of the radiative transfer equation for 
emergent intensity: 

/y = - {Sye~^Xvdz=- fcidz, (1) 

MJz MJz 

where is the source function, which is defined as the ratio 
between the emissivity to the opacity of the atmosphere; Ty is 
the monochromatic optical depth and the integrand C/ is the 
so-called intensity contribution function, which represents the 
intensity emanating from height z. 

4.1. Ila 6563 A line emission 

The Ha line is one of the most commonly observed lines 
that allow us to study the chromospheric response during 
a solar flare. It is sensitive to the fiux of non-thermal 
electrons precipitating to the chromosphere (jSyestkal [19761: 

iKasparovd et al .120091: iRubio da Costa et al.l2Q15h and a com¬ 

plex line formed in a broad height range, from the upper pho¬ 
tosphere to the lower chromosphere, making its interpretation 
non-trivial. 

The top panels of Figure [6] shows the temporal evolution 
of the resulting excess Ha line profile, where the quiet-Sun 
emission at t = 0 s has been subtracted to emphasize the 
changes during the fiare and the differences between Runs PL 
(red line) and SAl (blue line). In both runs, the line core 
presents a flattening due to the sudden behavior change of the 
source func tion in a very thin atmospher ic layer, as previously 
reported bv iRubio da Costa et al.l (l2015h . The line is broader 
and the intensity at the core is stronger for PL most of the 
time, while SAl shows stronger asymmetry (mostly blueshifts 
associated with upfiows) due to the higher velocities in the 
chromosphere. Close to t^ax = 60 s the profile presents a 
stronger blueshifted peak in the line center and redshifts in 
the wings, indicating plasma upfiows (evaporation) from the 
chromosphere to the transition region and downflows in the 
chromosphere. 

The above behavior can be better understood by examining 
the contribution function, as shown in Figure [71 We find that 
the photons in the wings of the line originate from the low 
chromosphere (at z = 0.12 Mm), where the plasma velocity is 
almost zero and almost constant in time. Therefore the wing 
emission changes little in time and is very similar in both runs. 
In contrast, the line core is formed in the transition region, 
whose height changes in time and where the plasma velocity 
suddenly increases due to upfiows. For SAl, this region is 
situated at a lower height and thus the overall height range of 
Ha line formation is narrower than that of PL. The PL model 
has lower velocities and covers a broader formation height 



0 20 40 60 80 100 

Time (seconds) 

Figure 8. Ha, He II 304 A, Ca II K and Si IV light curves. The emission at 
the initial time has been subtracted. The dashed line indicates tmax- [A color 
version of this hgure is available in the online journal.]. 

range at the line core. Moreover, the plasma temperature and 
density play an important role in the formation of the line and 
might ex plain the differences in the core emission between the 
two runs (iLeenaarts et akllMT^ . 

Integrating the intensity along the Ha line profile for a 
width of A/I =3.8 A yields the light curve shown in Fig¬ 
ure Oa). As can be seen, for SAl the Ha emission is directly 
correlated with the evolution of the fiux of electrons, peaking 
only 4 s after t^ax- However, for PL, the emission exhibits a 
plateau after 47 s, peaking exactly at tmax- 

4.2. He II304 A line emission 

The He II 304 A line is an optically thick line associ¬ 
ated with the transition '^S\i 2 - 2p Under ioniza¬ 
tion equilibrium conditions, in general, it is formed in a nar¬ 
row layer in the transition region, corresponding to tempera- 
tures of ^ 5 - 8 X 1 0"^ K (lO’Dwver et al.ll20T^ . However, 
iGolding et all (120141) note d that the det ails of its formation is 
still not well understood. IZirinI (Il988h showed that the pho¬ 
toionization and recombination processes play an important 
role in the He II 30 4 A emission. iJordanI (Il975h (see also 
lAndretta et al.l (l2003h ) found that the observed intensity of the 
helium-lines is higher than that estimated from observations 
of other EUV lines. They noted that the inclusion of high- 
energy electrons and cold ions in the n iodels enhances the 
synthe tic intensities. Based on this idea, iLaming & Feldman! 
(Il992h modeled the He I and He II emission during the im- 
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Figure 9. He II 304 A line profiles at r = 5, 20, 60 and 120 s for runs PL (red line) and SAI (blue thick line), similarly as Figure|6] 
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Figure 10. Intensity contribution function for He II 304 A as a function of frequency and height, similarly as Figure[7] Left panel: SAI Run; right panel: PL Run. 
The pink circle (in the SAI case - top) and triangle (in the PL case - bottom) mark the average formation height of the line center. The color scale of each panel 
is shown in its respective colorbar. Note that the line profiles are not in scale. [A color version of this figure is available in the online journal.]. 


pulsive phase of a flare by including the effects of high-energy 
electrons. 

In Figure [9] we note large differences between PL and SAI, 
indicating that the He II 304 A emission is very sensitive to 
the spectra of the non-thermal electrons. In fact, the intensity 
in PL changes by two orders of magnitude within 9 s of the 
simulation. Whenever there are large velocity variations in the 
formation region, the line profile shows strong asymmetries. 
In particular, the line profile of the PL model at f = 20 s shows 
a strong blueshift due to the evaporation of material to upper 
layers. 

Studying the monochromatic optical depth (green line in 


Figure (To]), we And that the photons of the He II 304 A line 
are emitted in a very narrow region situated at the bottom 
of the transition region. The height range of the formation 
region varies between 5 and 300 km for the PL model, be¬ 
ing even narrower (between 2 and 5 km after t = 10 s) for 
the SAI model. Figure [TT] shows the He III ionization frac¬ 
tion as a function of height and temperature. It indicates that 
the He II atoms become completely ionized within a narrow 
height range in the transition region. At r = 20 s, the PL model 
presents a gradual change of the ionization fraction over a 
broad height range, which is related to the broad formation 
heights shown in Figure [T0| at this time. 

From the intensity contribution function of Figure [T0| we 
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Figure 11. He III ionization fraction as a function of height (top) and tem¬ 
perature (bottom) for Runs SAl (left) and PL (right) at r = 5, 20, 60 and 
120 s. The pink circles (in the SAl case - left) and triangle (in the PL case - 
right) mark the height and temperature formation of the line center as shown 
in Figure[To] [A color version of this figure is available in the online journal.]. 

find that most of the intensity comes from +0.02 A from the 
line center, and thus the line center appears as a dip. The line 
profile is very sensitive to plasma velocity changes as shown, 
e.g. at f = 20 s for the PL model. The temperature range at the 
formation height is of the order of 1-2x10^ K (see the pink 
symbols in Figure [TT]). more than twice higher than the ex¬ 
pected values from lO’Dwver et all ( 201 Oh and consistent with 
the simulation results reported by iGolding et al.l (120141) . In 
addition. [Jordan et all (119931) found that the He II 304 A emis¬ 
sion is formed by collisional excitation in the quiet sun, but by 
the photoionization-recombination process in active regions 
and during fiares. This is consistent with the high tempera¬ 
ture values in the formation region revealed in our simulation. 

We note in passing that, as shown in Figure[8l the integrated 
intensity for the PL model experiences a rapid increase ait = 
11s and peaks at 23 s, which coincides with an increase in 
electron density at transition region temperatures. In contrast, 
such a peak is absent in the light curve of the SAl model, 
which exhibits more mild temporal variations. 


4.3. Ca IIK 3934 A line emission 

The Ca II K line is formed in the transition 2s -2p\, ab¬ 
sorbing photons at 3934 A. IViall (Il982b modeled the con¬ 
ditions at which the line profile is formed, assuming com¬ 
plete redistribution (CRD), requiring low electron densities 
(2 X 10^ ^ cm~^) and a relati vely low ionization degree of 
hydrogen. iPaletou et al.l (119931) improved the results by con¬ 
sidering the partial frequency redistribution approximation 
(PRD) in the model, getting higher intensity values. To model 
the Ca II K line profile it is important to keep in mind that 
the wings are formed under LTE con ditions, while the cor e 
is formed under non-LTE conditions (iLeenaarts et al.ll200^ . 
Another point to take into account in the modeling of this line 
is the formation height during the quiet Sun, which ranges 
from the chromosphere above the minimum temperature re¬ 
gion up to the transition region. 

The Ca II K line profile (Figure [12]) shows strong asymme¬ 


tries due to the large plasma velocities at these h eights, con¬ 
sistent with that observed by ICauzzi et al.l (119931) during the 
early phase of a solar flare. SAl presents stronger redshifted 
profiles (specially at f = 20 s) due to its higher velocities than 
PL, indicating that the Ca II K emission is sensitive to plasma 
velocities. 

The monochromatic optical depth (green line in Figure [13]) 
for Ca II K shows that the wings of the line are formed at 
z = 0.41 Mm in both runs, not changing in time. Instead the 
core is formed in the transition region, which is located at a 
lower height for the SAl model (e.g., see Figure[2]). As shown 
in Figure [13] the contribution function is weaker in the wings 
than in the line core, where the velocity changes rapidly with 
height. In fact, the plasma in the formation region is subject to 
downfiows and upfiows, which is manifested as the distortion 
in the line profile near the core. Nevertheless the wings still 
exhibit a symmetric shape. We also noted that the Ca IIK line 
wings are formed higher in the chromosphere than the Ha line 
for both runs. 

As Figure [Sfc) shows, the Ca IIK intensity in the PL model 
responds faster to the non-thermal electrons than in SAl, 
reaching higher intensity values at initial times; afterwards 
it follows a similar behavior as the fiux of the electrons, g- 
The fiux for the SAl model increases almost gradually, start¬ 
ing ai t = 12 s up to r = 65 s and decreases monotonically 
afterwards. 


4.4. Si IV1393 Aline emission 

The Si IV 1393 A line is an optically-thin line formed in 
a narrow region located in the upper chromosphere, which 
under ionization equilibrium conditi ons corresponds to a tem- 
perature range around (8 ±2) x 10"^ K. lDe Pontieu et al.l (120151) 
showed that the correlation between the non-thermal line 
broadening and the intensity is reproduced with simulations 
only when the non-equilibrium ionization is taken into ac¬ 
count. This is because the non-equilibrium ionization leads to 
the presence of Si^^ ions over a much wider rang e of temper¬ 
atures than under ionization equilibrium (see e.g. lOlluri et al.l 

[2on . 

We synthesized the intensity profile by using the CHI- 
ANTI abundances (iDere et al.l[l997l: iLandi et ^120131) and 
the plasma properties resulting from our atmosphere. The 
bottom panels of Figure [14] show the evolution of the Si IV 
line profile. This evolution can be described as an initial 
rapid blueshifted excursion due to chromospheric evaporation 
of material to the corona, followed by chromospheric com¬ 
pression causing redshifts and then blueshifts again in the late 
phase. 

Comparing the two runs, we note that, in general, PL 
presents larger blueshifts, while SAl has greater redshifts, es¬ 
pecially early in the simulation. This is due to their differ¬ 
ent plasma velocity distributions around the Si IV line forma¬ 
tion temperature shown in Figures [Ha) and (b). This trend 
is qualitatively co nsistent with the recent simulation result of 
[Testa et all ([20141) which shows that non-thermal electrons are 
required to produce blueshifted Si IV emission, while conduc¬ 
tive heating alone tends to produce only redshifted emission. 
The latter case resembles the strong heating in the corona 
rather than in the lower atmosphere in our Run SAl, due to the 
presence of its prominent quasi-thermal component of elec¬ 
trons. We also note differences in the large temporal varia¬ 
tions of the intensity between the two runs by up to orders of 
magnitude, the largest among all four emission lines studied 
here, as shown in Figures [14] and [S] Early in the simulation 
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Figure 12. Ca II K 3934 A line profiles at f = 5, 20, 60 and 120 s for runs PL (red line) and SAI (blue thick line), similarly as Figure[6| 
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Figure 13. Intensity contribution function for Ca II K as a function of frequency and height, similarly as Figure[7] Top panels: SAI Run; bottom panels: PL Run. 
All panels have the same color scale as indicated on the colorbar. [A color version of this figure is available in the online journal.]. 


during f =20 and 30 s, PL presents a broad peak in the in¬ 
tegrated intensity, more than 7 times higher than that of SAI 
(see Figure [ltd)). This is due to its larger electron density 
hump in the transition region, as shown in Figure [2l and thus 
greater emission measure. 

5. SUMMARY AND DISCUSSIONS 


The aim of this paper is to investigate the response of the 
solar atmosphere to the energy input by accelerated electrons. 
To achieve this, we have extended our earlier study (paper I) 
by including radiative transfer calculations. Specifically, we 
combined the Stanford Unified Acceleration-Transport code 
with the radiative hydrodynamics RADYN code. Our primary 
focus is to compare the results from the more realistic stochas- 























































12 



> 


CD 

> 


00 - 


CD 

cr 



- 0,1 0,0 0,1 - 0,1 0,0 0,1 - 0,1 0,0 0 , 1 ^ - 0,1 0,0 0,1 
Wavelength fronn line center [A] 

Figure 14. Si IV 1393 A line profiles at t = 5, 20, 60 and 120 s for runs PL (red line) and SAl (blue thick line), similarly as Figure|6] [A color version of this 
hgure is available in the online journal.]. 


tic acceleration model with that from an ad hoc power-law in¬ 
jection model and to obtain synthetic line emissions from our 
simulation results, which can provide new constraints for the 
particle acceleration mechanism. Our main findings are the 
following: 

1. In general, the temporal evolution of the atmosphere 
is determined not only by the energy fiux but also the 
spectral shape of non-thermal electrons. The results of 
our Runs PL (power-law injection) and SAl (stochastic 
acceleration) are in qualit ative agreement with those of 
Runs O and N of Paper I (iLiu et al.ll20Q9h . respectively. 
Stochastically accelerated electrons lead to stronger 
chromospheric evaporation with higher coronal temper¬ 
atures and plasma velocities (see Table [TJ because of 
their prominent quasi-thermal spectral component. Hy¬ 
drodynamic shocks form in both cases, lasting longer 
in PL, but appearing earlier in SAL 

2. The spatial distribution of the electron energy deposi¬ 
tion rate per unit volume is concentrated in the upper 
chromosphere for both runs (see bottom row of Fig¬ 
ure [2]). In terms of the electron heating rate per unit 
atmospheric mass, most of the directly deposited en¬ 
ergy is radiated away in the chromosphere for Run PL 
(Figure n~5]). However, due to the quasi-thermal spec¬ 
tral component. Run SAl has a significant amount of 
energy per unit mass directly deposited in the corona, 
rather than in the chromosphere (Figure [16]). Down¬ 
ward thermal conduction then carries the bulk of this 
energy to heat the transition region, where radiative 
loss is not as efficient as in the chromosphere. This 
is why chromospheric evaporation is much stronger in 
SAl than in PL, despite their identical total electron en¬ 
ergy fiux. 

3. In both cases, most of the energy exchange occurs in 
the lower atmosphere, and most of the net energy gain 
of the plasma is in the form of thermal energy, which 
dominates over ionization energy expenses. The ther¬ 
mal energy increase is about two orders of magnitude 
higher for SAl than for PL. 

4. For different SA runs of different peak electron en¬ 

ergy fiux the maximum upfiow velocity increases 
almost linearly with the while the maximum 


downflow velocity increases at a lower rate (see Ap- 
pendixjB]). The temperature peaks at t^ax is independent 

of g max- 

5. Surprisingly, explosive chromospheric evaporation 
with upfiow speeds of hundreds of km s“^ is present 
in our three SA runs with wide-ranging electron en¬ 
ergy fluxes ^max from 10^ to 10^^ erg s“^ cm“^. In 
contrast, explosive evaporation only occurred when 
'^m. ax >3 X 1 0^^ er g s“^ cm“^ in early simulations 
bv iFisher et H ] (fT^ . who used a power-law injec¬ 
tion with a spectral index of 4 and low-energy cutoff 
of 20 keV. Such energy fiux was widely quoted as the 
sole quantity that determines gentle vs. explosive evap¬ 
oration. Our result demonstrates that the shape of the 
electron spectrum is at least equally, if not more, im¬ 
portant as the total energy fiux, for the reasons noted 
above. This is co nsistent with the co nclusions indepen¬ 
dently reached bv iReep et al.l (120151) and ?. 

In order to study how the different models affect to the 
energy deposition in the lower chromosphere, we obtained 
the emission in several wavelength ranges formed at different 
heights from the upper photosphere to the transition region. 
According to the line emission for both runs, we find that the 
emission in Ha and Ca II K (two lines that cover the whole 
chromosphere) seems to be more dependent on the electron 
fiux variation than emission in other wavelengths, following 
a similar triangular temporal profile as the fiux of electrons, 
g. The characteristics of the specific emission lines are as 
follows: 

1. The Ha line is a broad line covering the chromosphere, 
from 0.12 Mm up to the transition region, where the 
temperature changes drastically. The SAl Run has 
higher chromospheric velocities and therefore stronger 
line-profile asymmetries than the PL Run. The core for¬ 
mation is located at a height where the plasma velocity 
changes drastically (from almost zero in the lower chro¬ 
mosphere to more than 50 km s“^), covering a broader 
region for the PL Run. Close to t^ax the profile presents 
a stronger blueshifted peak in the line center and red- 
shifts in the wings, indicating plasma upfiows (evapo¬ 
ration) from the chromosphere to the transition region 
and downflows in the chromosphere. The line intensity 
is temporally correlated with the electron energy fiux, 
peaking almost at its apex tmax- 
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2. Our synthetic He II 304 A line is formed within a tem¬ 
perature range from 1 x 10"^ to 2.5 x 10^ K, which 
is broader than previously reported for the quiet Sun, 
when this line is formed by collisional exci tation. How¬ 
ever, i n agreement with the results from [Jordan et aH 
(Il993h we find that, for active regions and during flares, 
the photoionization-recombination processes are more 
important. Thus the high temperatures at the formation 
height in our simulation are not unexpected. This line is 
formed at the bottom of the transition region, covering 
a height range between 5 and 300 km for the PL model, 
being even narrower for the SAl model. The shape 
of the line profile is strongly affected by the plasma 
velocity at the formation height, showing asymmetries 
at locations where large velocity gradients are present. 
The integrated intensity is strongly correlated with the 
electron density, as manifested in its large enhancement 
from 11 to 36 s for the PL model (see Figure [8tb)). 

3. The Ca II K line is formed in the chromosphere above 
the temperature minimum, at 0.41 Mm up to the transi¬ 
tion region, where the temperature changes drastically. 
The line profile shows strong asymmetries , consistent 
with the observations of ICauzzi et al.l (Il993h during the 
early phase of a solar flare. The Ca II K flux variation 
responds faster to the non-thermal electron heating in 
PL than in SAL 

4. The Si IV line is optically thin, formed in the up¬ 
per chromosphere. Run PL presents stronger chromo¬ 
spheric evaporation, greater up- and downflow veloc¬ 
ities and higher electron density values, resulting in 
seven times larger Si IV intensity than run SAL 

Our results have demonstrated that the emission line pro¬ 
files and light curves of the integrated intensities can provide 
useful plasma diagnostics as well as constraints for particle 
acceleration models. For example, the PL run exhibits similar 
behaviors of large increases in the light curves of both transi¬ 
tion region lines. He II 304 A and Si IV 1393 A, early during 
the simulation (see Figure [8]). Such large increases are in con¬ 
trast with the rather slow rise and even drop at certain times 
in the SAl case. Such predicted distinction in the temporal 
evolution can be checked against observations to shed light 
on the shape of non-thermal electrons. 

Our results also have important implications fo r sunquakes 
during solar flares (iKosovichev & Zharkovall9^ . Such seis¬ 


mic signals are believed to be caused by the response of 
the lower atmosphere to impulsive heating b y high-energy 
particles (iKosovichevI 120071: IZharkoval l2008h. among other 
proposed mechani sms (iLindsev & Doneall2008l: iDone 
iFisher et ^I2012^ ■ Specifically, sunquakes can be produced 
by hydrodynamic shocks propagating downward from the 
chromosphere to the photosphere at onsets of flares. In our 
simulations, we find that the hydrodynamic response to col¬ 
lisional heating results in initial downflows, especially in the 
SAl model, at velocities up to the order of 45 km s“^ within 
the first 10 s of the flare (see Figure [3tb)). Such high speed 
downflows are sufficient to create hydrodynamic shocks in 
the lower chromosphere and can be responsible for the sun- 
q uakes in the photo sphere, as demonstrated in the simulations 
of lZharko^ (120081) . 

The approach of combined particle acceleration-transport 
and radiative hydrodynamic simulation presented in this pa¬ 
per has opened a new chapter in modeling solar flare dy¬ 
namics. One important future improvement, among oth¬ 
ers, is to use nonparametric inver sion of acceleration pa¬ 
rameters from HXR ob servations (iPetrosian & Chenl 120101: 
IChen & Petrosianir2013l) as inputs to our model. This will al¬ 
low us to more accurately and realistically simulate the atmo¬ 
spheric evolution as well as the chromospheric and transition 
region emission lines, which can be checked against high- 
resolution spectroscopic observations, such as thos e from 
IRIS and DST/IBIS (e.g. iRuhio da Costa et alJl2()T^ . Such 
comparative investigations will provide promising new con¬ 
straints to the coupled processes in solar flares, including ra¬ 
diative transfer, hydrodynamics, and eventually particle accel¬ 
eration. 
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APPENDIX 

A. DISTRIBUTION OF ENERGY TERMS 

In this section, we examine the spatial distribution and temporal evolution of different terms of the energy equation (lAlb . in 
order to better understand their contribution to the internal energy of the plasma. 

I ^ I" (P + + Fj-detail) “ Q'mii ~ Qe ~ SxEUV — 0? (Al) 

Ot oz oz oz 

Here z is distance, t time, p density, v velocity, p pressure, e internal energy per unit mass, qy the viscous stress, and Fy 
and Fr are conductive and radiative fluxes. The specific energy terms include: viscous heating -qy{dvldz), compression work 
-p(dv/dz), conductive heating -dFddz, the heating function to create the initial atmosphere, collisional heating Qy by 
non-thermal electrons. X-ray plus EUV heating gxEUV, the thin radiative heating -dF^^thi^^ldz of elements treated in LTE and the 
radiative flux detail/calculated from solving the radiative transfer equation. Here Er detail is integrated over wavelength and 
includes about 34000 wavelengths in the temperature range between 0.1 and 10 MK, covering the continua from 1 to 40000 A 
and spectral lines of Hydrogen, Helium, Calcium and Magnesium, treated in Non-LTE. 
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Figure 15. Variations with height of different energy terms for PL at r = 5, 20, 60 and 120 s. Note that the X-axis is in logarithmic scale. Top: rate of change 
of major energy terms, including the total internal energy (black thick line), radiative heating (pink), conductive heating (cyan), electron heating (green), and 
compression work (red). Bottom: rate of change of the total internal energy (black), divided into thermal energy (red), hydrogen ionization and excitation energy 
(blue) and ionization and excitation energy for elements other than hydrogen (green). The insets show enlarged views near the transition region. [A color version 
of this figure is available in the online journal.]. 

Since the atmospheres evolves differently in the two runs, we discuss them individually. In order to avoid complex graphics, 
here we show only the energy terms with significant contributions to the total internal energy. 


AL PL Run 

The top panels of Figure [15] show the rate of change of different energies per unit mass (cf., heating rate per unit volume in the 
fourth row of Figure [2]): internal energy (black thick line), radiative heating (pink), conductive heating (cyan), electron heating 
(green), and compression work (red). Most of the energy exchange occurs in the lower atmosphere and the main contribution to 
the internal energy change is the electron heating rate. At f = 5 s and z = 1.55 Mm (within the chromosphere), the conductive 
and electron heating rates are the main source of heating, which are mostly balanced by the radiative loss term, explaining the 
small temperature changes at this time (see Figure [2]). 

At f = 20 s, deep in the chromosphere, electron heating Qe is largely balanced by radiative losses. However, in the transition 
region (near z = 2 Mm) and corona, it dominates over the loss terms and results in an increase in the internal energy and 
temperature. Around the chromospheric evaporation front at z = 3.7 Mm, the pressure work due to the negative velocity gradient 
leads to a local temperature enhancement (see Figure|2]). Thus, compensating conductive cooling. 

At later times, conductive heating in the transition region becomes important. For example, at f = 60 s and z = 1.07 Mm, the 
conductive heating rate increases in a height range of less than 150 m from almost zero to 3 x 10^^ (erg s“^ g“^). This is also the 
region where most of the non-thermal electron energy is deposited (see bottom row of Figure [2]). 

The bottom panel of Figure [15] shows the change rate of the total internal energy (black thick line) divided into the thermal 
energy of the material (red), the sum of ionization energy of hydrogen (blue) and the rest of the constituent atoms (green). As 
can be seen, the thermal energy change dominates at all times, but the total ionization energy of atoms other than hydrogen is 
important at early times. 


A2. SAl Run 

In the SAl model, the electron heating rate per unit mass, the major contributor to the total internal energy increase (see 
Figure [16]), is primarily located in the corona, rather than the transition region or chromosphere as in the case of the PL model. 
This results in the significant increase of the coronal temperature (see Figure [2]). At ^ = 5 s the pressure work balances the 
conductive cooling at z = 1.86 Mm (above the transition region), where the electron density is higher. 

At r = 20 s the bulk of pressure work has moved to the upper corona and is still balanced by conductive cooling. At z = 
1.125 Mm, where the transition region is located, the conductive heating rate increases sharply up to 4.4 x 10^^ (erg s“^ g“^) and 
is partly radiated away and balanced by the pressure work within a very narrow region, where the non-thermal electron energy 
is mostly deposited (see bottom row of Figure [2]). The electron heating in the lower atmosphere is negligible in comparison with 
the other terms and it is mostly radiated away. 

As the atmosphere evolves, most of the conductive energy in the transition region is radiated away and very little is left to heat 
the plasma. This energy exchange occurs in a very narrow layer, less than 10 m thick, at 0.866 Mm, where the conductive 
heating rate increases from almost 0 to 5 x 10^^ erg s“^ g“^ (see the inset in Figure [T6tc)). The atmosphere is already ionized 
in this region and does not radiate efficiently anymore, thus increasing the total internal energy and temperature. Note that for a 
better comparison, the last two columns of Figure [16] do not show the total length of the peak. 
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Figure 16. Same as Figure[T5] but for Run SAL Note that the X-axis is in logarithmic scale and that at r = 60 and 120 s, the peaks at the transition region height 
are truncated in the main panels but are displayed in full in the insets. [A color version of this hgure is available in the online journal.]. 


As shown in the bottom panels, most of the total energy gain (black thick line) is used to increase the thermal energy (red line), 
as in the case of PL. However, these energy change rates are one to two orders of magnitude higher than those in PL. 

B. EFFECTS OE THE AMPLITUDE OE THE ELECTRON ELUX ^max 

Here we present the result of SA runs of different peak energy fluxes of the non-thermal electrons, gmax. but with the same 
triangular-shaped temporal profile of g(0 as in Runs SAl and PL. The characteristics of these runs are shown in Table [TJ As 
expected, a higher electron flux leads to higher maximum upflow and downflow speeds and temperatures which occur at earlier 
times. The only exception is the maximum temperature that almost always occurs at the time of Smax, be., f = 60 s. A factor of 
20 increase in gmax from 5.7 x 10^ erg s“^ cm“^ (Run SA3) to 1.2 x 10^^ erg s“^ cm“^ (Run SAl) only results in a factor of two 
increase in the maximum upflow speed Vmax from 377 to 750 km s“\ a factor of 50% increase in the maximum downflow speed 
Vmin from -28 to -41 kms“^ and a factor of three increase in the maximum temperature from 8 to 23 MK. These results are in 
line with those of lLiul[2008l see his Table 7.1 for simulations using the same spectrum of stochastically accelerated electrons as 
we have adopted here. 


t = 5 secs t = 20 secs t = 60 secs t=120 secs 
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Figure 17. Spatial distribution of the electron heating rate Qe{z) for Runs SAl, SA2 and SA3 (see Table[TJ at selected times, r = 5, 20, 60 and 120 s. [A color 
version of this figure is available in the online journal.]. 


iFisher et ^ (Il985h found that an electron energy flux of < 10^^ erg s“^ cm“^ can produce so-called “gentle” chromospheric 
evaporation with upflow speeds of tens of km s“\ while an energy flux of > 3 x 10^^ erg s“^ cm“^ would produce “explosive” 
evaporation with upflow speeds of hundreds of km s“^ and downflow speeds of tens of km s“^. As shown in Table [TJ all the SA 
runs have energy flux below this threshold and yet all exhibit characteristics of “explosive” evaporation. This demonstrates that 
the energy flux is not the only factor that dictates the atmospheric response to electron heating; the electron spectral shape is 
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also critical because it determines the spatial distribution of the energy deposition. The key difference is that iFisher et al.l (Il985h 
injected a single power-law electron spectrum of an index of d = 4 and cutoff energy of 20 keV, which is comparable to that 
in our PL run. Such non-thermal electrons directly deposit most of their energy deep in the chromosphere, which is largely 
radiated away (see Figure [15]). In our SA runs, the electron spectrum has a prominent quasi-thermal component and thus results 
in heating primarily in the corona (in terms of heating rate per unit mass, see Figure [16]). Then thermal conduction becomes 
the primary heating agent in the lower atmosphere, mainly in the transition region where radiative loss is not as strong as in the 
chromosphere. As a result, there is relatively more energy left to heat the plasma and drive chromospheric evaporation. This is 
similar to the case in Paper I, although radiative transfer was not included there. These result indicates that the stochastically 
accelerated electrons, because of their preferential heating in the corona rather than directly in the chromosphere, are more 
efficient in driving chromospheric evaporation than power-law electrons with a cutoff energy. This also indicates that, in addition 
to direct collisional heating by non-thermal electrons as commonly believed, thermal conduction can play an iin portant role in 
driving chromospheric evaporation, as observed in some flares fe.g.. IZarro & Lem^ll988l: [Battaglia et al.ll200^ . This is even 
more so for weak flares, such our SA3 run. 

In the case of faint flares (Run SA3) most of the electron heating rate, Qe, is transformed into conductive heating and radiated 
away along the atmosphere, explaining why the electron heating in Figure [T7] (green line) remains almost constant in the corona. 
By increasing the electron flux g, in Runs SAl and SA2 (blue and yellow lines), the heating increases at the top of the loop with 
time. 


Bl. Lines emission 
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Figure 18. Light curves of the Her 6563 A, He II 304 A, Ca II K 3934 A, and Si IV 1393 A lines for Runs SAI-SA3. The emission at the initial time has been 
subtracted. The vertical dashed line indicates [A color version of this hgure is available in the online journal.]. 


Figure [18] shows the light curve of the Ha, He II 304 A, Ca II K and Si IV 1393 A lines for the three SA models, where 
the emission of the quiet Sun (at the initial time) has been removed. In general, for higher peak electron energy fluxes, the line 
intensity increases faster with time, because of more rapid hydrodynamic response of the atmosphere, and reaches higher peak 
values. This is the case for the Ha and Ca II K fluxes, which are positively correlated with the electron flux, peaking around tmax- 
However, the He II 304 A intensity in Run SAl decreases in the middle portion of the flare duration, showing anti-correlation 
in time with the electron flux. Each of the He II light curves of the three runs also shows distinct behaviors. In Table [T] we 
can see that the maximum temperature decreases linearly with the flux of the injected electrons. Therefore the recombination- 
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photoionization processes may play a less important role in the formation of the line, giving preference to the collision excitation 
processes. This may explain why the shape of the He II 304 A line profiles and the flux changes with g- 
The temporal evolution of the Si IV 1393 A flux shows an earl ier b ump associated with the increase of the electron density in 
the chromosphere and transition region. As discussed in Section [441 an increase of the electron density leads to an increase of 
the emission measure at this wavelength range and thus the increase in intensity of this optically-thin line. 

REFERENCES 


Abbett, W. R & Hawley, S. L. 1999, ApJ, 521, 906 

Allred, J., Kowalski, A., & Carlsson, M. 2015, in AAS/AGU Triennial Earth-Sun Summit, Vol. 1, AAS/AGU Triennial Earth-Sun Summit, 30207 

Allred, J. C., Hawley, S. L., Abbett, W. R, & Carlsson, M. 2005, ApJ, 630, 573 

Andretta, V., Del Zanna, G., & Jordan, S. D. 2003, A&A, 400, 737 

Battaglia, M., Fletcher, L., & Benz, A. O. 2009, A&A, 498, 891 

Carlsson, M. & Stein, R. F. 1992, ApJ, 397, L59 

—. 1997, ApJ, 481,500 

Cauzzi, G., Falchi, A., Falciani, R., & Smaldone, L. A. 1993, Advances in Space Research, 13, 311 
Chen, Q. & Retrosian, V. 2013, ApJ, 777, 33 

De Rontieu, B., McIntosh, S., Martinez-Sykora, J., Reter, H., & Rereira, T. M. D. 2015, ApJ, 799, L12 

De Rontieu, B., Title, A. M., Lemen, J. R., Kushner, G. D., Akin, D. J., Allard, B., Berger, T., Boerner, R, Cheung, M., Chou, C., Drake, J. E, Duncan, D. W., 
Freeland, S., Heyman, G. E, Hoffman, C., Hurlburt, N. E., Lindgren, R. W., Mathur, D., Rehse, R., Sabolish, D., Seguin, R., Schrijver, C. J., Tarbell, T. D., 
Wiilser, J.-R, Wolfson, C. J., Yanari, C., Mudge, J., Nguyen-Rhuc, N., Timmons, R., van Bezooijen, R., Weingrod, L, Brookner, R., Butcher, G., Dougherty, 
B., Eder, J., Knagenhjelm, V., Larsen, S., Mansir, D., Rhan, L., Boyle, R, Cheimets, R N., DeLuca, E. E., Golub, L., Gates, R., Hertz, E., McKillop, S., Rark, 
S., Rerry, T., Rodgorski, W. A., Reeves, K., Saar, S., Testa, R, Tian, H., Weber, M., Dunn, C., Eccles, S., Jaeggli, S. A., Kankelborg, C. C., Mashburn, K., 
Rust, N., Springer, L., Carvalho, R., Kleint, L., Marmie, J., Mazmanian, E., Rereira, T. M. D., Sawyer, S., Strong, J., Worden, S. R, Carlsson, M., Hansteen, 
V. H., Leenaarts, J., Wiesmann, M., Aloise, J., Chu, K.-C., Bush, R. L, Scherrer, R H., Brekke, R, Martinez-Sykora, J., Lites, B. W., McIntosh, S. W., 
Uitenbroek, H., Okamoto, T. J., Gummin, M. A., Auker, G., Jerram, R, Rool, R, & Waltham, N. 2014, Sol. Rhys., 289, 2733 
Dere, K. R, Landi, E., Mason, H. E., Monsignori Fossi, B. C., & Young, R. R. 1997, A&AS, 125, 149 
Donea, A. 2011, Space Sci. Rev., 158, 451 

Dorfi, E. A. & Drury, L. O. 1987, Journal of Computational Rhysics, 69, 175 
Emslie, A. G. 1978, ApJ, 224, 241 
—. 1981, ApJ, 249,817 

Fisher, G. H., Bercik, D. J., Welsch, B. T., & Hudson, H. S. 2012, Sol. Rhys., 277, 59 
Fisher, G. H., Canfield, R. C., & McClymont, A. N. 1985, ApJ, 289, 414 
Golding, T. R, Carlsson, M., & Leenaarts, J. 2014, ApJ, 784, 30 

Gustafsson, B. 1973, A FORTRAN Rrogram for Calculating ’Continuous’ Absorption Coefficients of Stellar Atmospheres (Uppsala: Lalidstingets Verkstader) 
Hamilton, R. J. & Retrosian, V. 1992, ApJ, 398, 350 
Holman, G. D. 1985, ApJ, 293, 584 

Hoyng, R, Duijveman, A., Machado, M. E., Rust, D. M., Svestka, Z., Boelee, A., de Jager, C., Frost, K. T., Lafleur, H., Simnett, G. M., van Beek, H. E, & 
Woodgate, B. E. 1981, ApJ, 246, L155 
Johannesson, A., Marquette, W., & Zirin, H. 1995, Sol. Rhys., 161, 201 
Johannesson, A., Marquette, W. H., & Zirin, H. 1998, Sol. Rhys., 177, 265 
Jordan, C. 1975, MNRAS, 170, 429 

Jordan, S. D., Thompson, W. T., Thomas, R. J., & Neupert, W. M. 1993, ApJ, 406, 346 

Kasparova, J., Varady, M., Heinzel, R, Karlicky, M., & Moravec, Z. 2009, A&A, 499, 923 

Kosovichev, A. G. 2007, ApJ, 670, L65 

Kosovichev, A. G. & Zharkova, V. V. 1998, Nature, 393, 317 

Kundu, M. R., White, S. M., Gopalswamy, N., & Lim, J. 1994, ApJS, 90, 599 

Laming, J. M. & Feldman, U. 1992, ApJ, 386, 364 

Landi, E., Young, R. R., Dere, K. R, Del Zanna, G., & Mason, H. E. 2013, ApJ, 763, 86 
Leach, J. & Retrosian, V. 1981, ApJ, 251, 781 

Leenaarts, J., Carlsson, M., & Rouppe van der Voort, L. 2012, ApJ, 749, 136 

Leenaarts, J., Rutten, R. J., Stitterlin, R, Carlsson, M., & Uitenbroek, H. 2006, A&A, 449, 1209 

Lin, R. R. 1985, Sol. Rhys., 100, 537 

Lindsey, C. & Donea, A.-C. 2008, Sol. Rhys., 251, 627 

Liu, S., Retrosian, V., & Mason, G. M. 2004, ApJ, 613, L81 

—. 2006, ApJ, 636, 462 

Liu, W. 2008, Solar Flares as Natural Rarticle Accelerators: A High-energy View from X-ray Observations and Theoretical Models (VDM Verlag Dr) 

Liu, W., Chen, Q., & Retrosian, V. 2013, ApJ, 767, 168 
Liu, W., Retrosian, V., & Mariska, J. T. 2009, ApJ, 702, 1553 
Mariska, J. T., Emslie, A. G., & Li, R. 1989, ApJ, 341, 1067 
McTiernan, J. M. & Retrosian, V. 1990, ApJ, 359, 524 
Milkey, R. W. & Mihalas, D. 1973, ApJ, 185, 709 
Miller, J. A. 1997, ApJ, 491, 939 
Milligan, R. O. & Dennis, B. R. 2009, ApJ, 699, 968 

O’Dwyer, B., Del Zanna, G., Mason, H. E., Weber, M. A., & Tripathi, D. 2010, A&A, 521, A21 

Olluri, K., Gudiksen, B. V., & Hansteen, V. H. 2013, ApJ, 767, 43 

Raletou, E, Vial, J.-C., & Auer, L. H. 1993, A&A, 274, 571 

Rark, B. T. & Retrosian, V. 1995, ApJ, 446, 699 

Retrosian, V. 2012, Space Sci. Rev., 173, 535 

Retrosian, V. & Chen, Q. 2010, ApJ, 712, L131 

Retrosian, V. & Donaghy, T. Q. 1999, ApJ, 527, 945 

Retrosian, V., Donaghy, T. Q., & McTiernan, J. M. 2002, ApJ, 569, 459 

Retrosian, V. & East, W. E. 2008, ApJ, 682, 175 

Retrosian, V. & Liu, S. 2004, ApJ, 610, 550 


18 


Potzi, W., Temmer, M., Veronig, A., Hirtenfellner-Polanec, W., & Baumgartner, D. 2013, in EGU General Assembly Conference Abstracts, Vol. 15, EGU 
General Assembly Conference Abstracts, 1459 

Ramaty, R. 1979, in American Institute of Physics Conference Series, Vol. 56, Particle Acceleration Mechanisms in Astrophysics, ed. J. Arons, C. McKee, & 

C. Max, 135-154 

Reep, J., Bradshaw, S., & Alexander, D. 2015, ArXiv e-prints 

Rubio da Costa, E, Kleint, L., Petrosian, V., Sainz Dalda, A., & Liu, W. 2015, ApJ, 804, 56 
Sakao, T. 1994, PhD thesis, (University of Tokyo), (1994) 

Sironi, L. & Spitkovsky, A. 2009, ApJ, 707, L92 

Testa, R, De Pontieu, B., Allred, J., Carlsson, M., Reale, K, Daw, A., Hansteen, V, Martinez-Sykora, J., Liu, W., DeLuca, E. E., Golub, L., McKillop, S., 
Reeves, K., Saar, S., Tian, H., Lemen, J., Title, A., Boerner, R, Hurlburt, N., Tarbell, T. D., Wuelser, J. R, Kleint, L., Kankelborg, C., & Jaeggli, S. 2014, 
Science, 346, B315 

Tsuneta, S. & Naito, T. 1998, ApJ, 495, L67 

Svestka, Z. 1976, Several Solar Aspects of Flare-Associated Particle Events, ed. P. S. McIntosh & M. Dryer, 141-162 
Vial, J. C. 1982, ApJ, 253, 330 

Woods, T. N., Eparvier, F. G., Hock, R., Jones, A. R., Woodraska, D., Judge, D., Didkovsky, L., Lean, J., Mariska, J., Warren, H., McMullin, D., Chamberlin, R, 
Berthiaume, G., Bailey, S., Fuller-Rowell, T., Sojka, J., Tobiska, W. K., & Viereck, R. 2012, Sol. Phys., 275, 115 
Zarro, D. M. & Lemen, J. R. 1988, ApJ, 329, 456 
Zharkova, V. V. 2008, Sol. Phys., 251, 641 
Zirin, H. 1988, Astrophysics of the sun 


